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Abstract. This paper presents and discusses the application of blind 
source separation to astrophysical data obtained with the WMAP satel- 
lite. Blind separation permits to identify and isolate a component com- 
patible with CMB emission, and to measure both its spatial power spec- 
trum and its emission law. Both are found to be compatible with the 
present concordance cosmological model. This application demonstrates 
the usefulness of blind ICA for cosmological applications. 



1 Introduction 

The detection of Cosmic Microwave Background (CMB) fluctuations, and the 
measurement of their spatial power spectrum, has been over the past three 
decades subject of intense activity in the cosmology community. 

The CMB, discovered in 1965 by Penzias and Wilson, is a rehc radiation 
emitted some 13 billion years ago, when the Universe was about 370.000 years 
old. This radiation emits, as a function of wavelength, as an almost perfect black- 
body at a temperature of 2.726 Kelvin. Small fluctuations of this temperature 
tracing the seeds of the primordial homogeneities which gave rise to present large 
scale structures as galaxies and clusters of galaxies, are expected to exist, at the 
level of a few tens of micro Kelvin. 

The importance of measuring these fluctuations to constrain cosmological 
scenarios describing the history and properties of our Universe is now well es- 
tablished. The measurement of their statistical properties (and in particular 
of their spatial power spectrum) permits to drastically constrain the cosmo- 
logical parameters describing the matter content, the geometry, and the evolu- 
tion of our Universe |Jungman et al. 1996| . The accuracy required for precision 
tests of the cosmological scenarios, however, is such that it is necessary to dis- 
entangle in the data the contribution of several distinct astrophysical sources, 
all of which emit radiation in the frequency range used for CMB observations 
[Bouchet fc Gispert 1999| . In addition, all measurements being noisy, it is nec- 
essary to identify, characterize, and remove as well as possible contributions due 
to instrumental noise. 

CMB emission peaks in the millimeter-wave domain of the electromagnetic 
spectrum. In this wavelength range however, several astrophysical sources of 



radiation contribute to the total emission. Some of them arc quite well known 
(or modelled), and others not so well. It is of utmost importance then, for the 
interpretation of the observations, to separate the various emissions. This can 
be achieved (to some extent) by the joint processing of various observations 
performed at various wavelengths. 

To first order, the observation at a given wavelength A can be modeled as a 
linear superposition of a number of sources. A sequence {X(x, y)} of observations 
— J maps obtained at several wavelengths — is modeled as an instantaneous 
mixture of K independent sources contaminated by independent noise N(a;,y): 

X{x,y) = AS(a;,y) +N(x,y) 

where A is an unknown J x K matrix. 

As a strong prediction of the standard big bang model of cosmology, CMB 
anisotropies are expected to have a spectral emission law given, to first order, 
by the derivative (with respect to temperature) of the blackbody law, taken at 
the temperature of the CMB: 
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where B^(T) is the Planck law for the emission of a blackbody at temperature T. 
Blind component separation, which permits to estimate the emission law of sky 
components, provides a unique tool for checking this prediction. In the following, 
we apply a blind component separation method on WMAP data. 



2 WMAP data 
2.1 WMAP 

The WMAP space probe, launched by NASA in 2001, is a large telescope 
for imaging the total emission of the sky at 5 different wavelengths (or fre- 
quency channels), with a resolution ranging from about 0.2 degrees to 0.9 de- 
grees (limited by diffraction), and with full sky coverage {An of the sphere) 
[ Bennett et al. 2003j . 

Data taken with the WMAP instrument has been made available to the 
scientific community after one year of proprietary period, and is freely available 
on a dedicated web site.^ In its most usable format, the data consists of a set 
of 10 maps obtained by different detector pairs: four maps at 3.2 mm (94 GHz), 
two at 4.9 (61 GHz), and one at 7.3, at 9.1, and at 13.0 mm (41, 33 and 23 
GHz). The data is provided in the Healpix pixellisation format of the sphere.^ 
Data from the various detectors observing the sky at the same frequency can be 
combined to yield a single map of the sky at that frequency. In this analysis, we 
use the data from individual detector pairs. 

* http://lambda.gsfc.nasa.gov/ 

^ http: / /www. eso.org/ science /healpix/ 



2.2 Astrophysical emissions 

In the millimeter waveband, the total sky emission can be considered as the su- 
perposition of two main classes of sources: diffuse emissions and compact sources. 

Diffuse emissions are due to extended processes on the sky. For all of these 
processes, the observed photons are emitted by large scale distributions of emit- 
ters - gases of particles. Among these, one expects emission from the hot plasma 
in the early Universe (the CMB), greybody emission from cold dust particles in 
the galaxy (dust emission), emission from relativistic electrons spiraling in the 
galactic magnetic field (synchrotron), free-free (Bremsstrahlung) emission from 
ionized gaz. Compact source emission is due to distant objects (galaxies, quasars, 
clusters of galaxies) , the angular size of which is smaller than the resolution of 
the instrument. 

For each physical source of emission, the total observable intensity /(A, 0, 0), 
depends on the physical properties of the source (dumpiness, emission process, 
temperature, etc.). Here, A is the wavelength and (6*, </>) is the direction on 
the celestial sphere. For most diffuse emissions, we can approximate /(A, 9, (j)) ~ 
A(A)S(6', (/)), the approximation being expected to be excellent (less than 0.1 %) 
for some components (e.g. CMB) and not so good (at the 10% level) for others 
as dust emission. 

2.3 Spherical harmonics 

It is convenient and usual to decompose the sky brightness (traditionnally ex- 
pressed as AT {9, if) in units of temperature fluctuation in micro Kelvin), as an 
expansion: 

z\r(0,^) -^a,^,„r,,™(0,¥>) (2) 
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where the functions y£_m(0, i^) are the usual spherical harmonics. For a stationary 
isotropic random Gaussian field on the sphere, the coefficients ai^m on the sphere 
are such that ^{ai^mO-i' ,m') = Cidwdmm' ■ On the spherical harmonics basis, a 
noisy linear mixture of components writes 

X(^, m) = AS(^, m) -(- N(£, m) (3) 

The WMAP observations used in this application are decomposed on spher- 
ical harmonics, and the component separation method is applied on the set 
{X(€, m)}, a 10 X 1 vector of coefficients of the spherical harmonics expansion 
of the WMAP observations. 

3 Method 

3.1 Spectral Matching 

Following the method described in [Delabrouille et al. 2003| , we compute average 
values of the covariance of the observations: 

Rx(g)-- V X(€,m)X(^,m)t {q = l,...,Q) (4) 

72 



where q labels a domain Vq of {£, m) space and Uq = card(Pq). For our applica- 
tion, it is customary to estimate band averages of the power spectra in domains 
£i < £ < £2, with m taking all available values for each £, i.e. —£< m < £. 

If the observations X match the multi-component model of equation |3| then 
Rx(9) is an estimator of Rx(9), the latter being itself given by 

Rx(g) = ARs(g)At+RN(g) (<7=1,...,Q) (5) 

with Rs(9) and Rn(9) defined similarly to Rx('i')- The parameters of interest 
(mixing matrix and band-averaged spatial power spectra of components) are 
precisely A and Rs(g), the latter being diagonal for uncorrelated components. 
Noise spatial power spectra, being of no particular use for astrophysics, are 
nuisance parameters which must be estimated, so that the full set of parameters 
of the data is 6 = {A, Rs, Rn}- 

The estimation of these parameters of the model is done by maximizing an 
approximation of the likelihood (in the Whittle approximation), which (up to 
irrelevant factors) boils down to minimizing 

Q 

</.(0) = ^n, Z?(Rx(g),Rx(<z;e)) (6) 

9=1 

D{-, ■), the KuUback-Leibler divergence, is a measure of divergence between two 
positive n X n matrices, defined by 

i:'(Ri,R2) = ^ [tr (RiRs"^) - logdet(RiR^^) - n] (7) 

Optimization is achieved with a dedicated algorithms based on the EM algo- 
rithm, followed by few steps of a descent method (BFGS algorithm) to speed-up 
the convergence [Cardoso fc Pham] . 

3.2 Estimation of Component maps 

Component maps are estimated by applying to the observations X a Wiener 
filter W — [A*R^^A -I- Rg ^]~-^A*Rj;j^, where estimated values A, Rn and Rg 
of the parameters are used. The estimated components are obtained then by: 

§ = [A*Rn^A + Rs^]"^A*Rn1X (8) 

In the limiting case where noise is small as compared to component signals, 
Rg ^ is negligible. If in addition estimates of the mixing matrix is good (A ~ A) 
the application of the Wiener filter to the observations yields 

S~S-K[A*Rn1A]-1A*Rn^N (9) 

which is an unbiased estimate of S. In poor signal to noise regimes, the signal is 
attenuated to suppress noise contamination in the reconstructed maps. 



3.3 Dealing with real data 



The observations of astrophysical emissions are more complicated than the sim- 
ple model of equation O in several respects. 

— Known galactic emissions (dust, synchrotron, free- free) are positive and con- 
centrated towards the galactic plane; 

— In addition, galactic component emissions may be correlated for physical 
reasons (thus violating the independence assumption); 

— The emission law (mixing matrix elements) of galactic components is known 
to be slightly position dependent and/or scale dependent; 

— Point sources are not modeled satisfactorily by a linear mixture, and their 
distribution is poorly represented by the Whittle approximation; 

— Instrumental noise is not stationary over the sky, which reduces somewhat 
the efficiency of parameter estimation; 

— Detectors observe the sky with a detector dependent angular resolution. This 
effect can be modeled, to first order, as a detector-dependent multiplicative 
term Bi in harmonic space. 

For the present analysis, we deal with these complications in the simplest pos- 
sible way. We decide to focus the present analysis on checking whether WMAP 
data contains a component which has both a CMB power spectrum and the ex- 
pected emission law for temperature fluctuations of the cosmological blackbody. 
Hence, we cut out from the analysis those regions of the sky which are most 
contaminated by galactic emission. Also, we cut out from the analysis small 
regions around known bright point sources (see next subsection). This reduces 
drastically the impact of the first four items in the list. We ignore noise non sta- 
tionarity. Irregularities are smooth enough on the sky to induce only small (and 
small range) correlations in harmonic space. The impact of the multiplicative 
term Bg in harmonic space is corrected for beforehand. 

We restrict the analysis to multipole t values between 10 and 400, which 
encompasses all of the first acoustic peak observed by WMAP using the highest 
frequency channels only. The upper limit permits to minimize the impact of 
a background of residual unresolved point sources and of uncertainties in the 
modeling of the beam of the instrument. 

3.4 Mask 

We apply a mask to the observations which sets to regions of the sky con- 
taminated by strong emission which probably are not well modeled as a linear 
mixture of uncorrelated sources. These regions are essentially a set of small cuts 
centered on known bright point sources, and a large cut through the map to 
exclude the galactic plane region, highly contaminated by a complex set of cor- 
related foreground emissions. This cut permits to isolate a region (the major 
fraction of the sky) where the emission is thought to be dominated by CMB 
fluctuations. 



It is possible to perform the SMICA analysis using the decomposition in 
spherical harmonics of cut — sky data. The analysis, then, is performed on a 
fraction of the sky (at galactic latitude higher than 20 degrees, with point sources 
removed, and with small additional high-level galactic emission regions cut). 

On such incomplete sky, some technical issues have to be taken into account. 
First, the spherical harmonics do not form an orthonormal basis, but his has lit- 
tle impact on the estimation in large bands oi £. Second, because of the mask, the 
stationarity hypothesis which guarantees the decorrelation of spherical harmon- 
ics expansion coefficients of the components is not quite satisfied, which results 
in possible loss of efficiency. As a final comment on masking, we note that the 
spectra estimated on the cut sky have to be corrected to take into account the 
impact of this "window function" . We apply a correction using the Master for- 
malism of Hivon et al. | Hivon et al. 2002) to correct from partial sky coverage 
effects. 

4 Results 

The application of SMICA to WMAP data gives the following results: 

— Two components are clearly identified, one of which is CMB, and the second 
emission from the high latitude galaxy. 

— A third component is identified marginally, which may correspond to large 
scale variations of galactic emission law, 

— The identified CMB component has an emission law compatible with the 
expected derivative of a blackbody, to excellent accuracy. Its power spectrum 
displays a peak around £ = 200 (first acoustic peak) compatible with the 
measurement announced by the WMAP team. 

— The galactic component has an emission law compatible with synchrotron 
emission in all WMAP channels (proportional to z^~^'^) except those at 90 
GHz, for which a significant excess is seen. This excess is interpreted as 
due to galactic dust emission correlated with the synchrotron (either for 
physical reasons or - more likely ~ by happenstance because of concentration 
of galactic emissions at low galactic latitude). 

Error bars for all quantities are computed from the Hessian at the point of 
convergence, with a correction factor l/\/{,fsky for part-sky coverage. 

5 Conclusion 

The SMICA blind component separation has been applied successfully to WMAP 
data. The analysis is performed on £ varying from about 10 to 400, which encom- 
passes all of the first acoustic peak of the CMB. We find that all 10 WMAP maps 
comprise a common astrophysical component identified here as CMB. The emis- 
sion law of this first component is compatible with CMB fluctuations to within 
about one per cent, and its spatial power spectrum in the range 10 < ^ < 400 
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Fig. 1. Left : Measured CMB fluctuation mixing matrix, compared to theoret- 
ical expectations. Data points at the same frequency have been shghtly offset 
for readabihty. This figure shows that the first component identified bUndly in 
the WMAP data, and seen consistently in all WMAP channels, has a CMB fiuc- 
tuation emission law to within less than one per cent (error bars here exclude 
WMAP calibration errors). Right: Measured CMB spatial power spectrum, com- 
pared to the theoretical model matching best WMAP team measurement (solid 
black line). This figure shows that the first component identified blindly in the 
WMAP data, and seen consistently in all WMAP channels, has a spatial power 
spectrum consistent with what has been measured by WMAP using only the 
two highest frequency channels. 



compatible with the measurement obtained by the WMAP team. The second 
component has an emission law proportional to except in the W (90 GHz) 

channel for which a significant excess is observed. This excess is thought to be 
due to dust emission contribution to the component. These results demonstrate 
the usefulness of blind component separation for astrophysical data analysis. 
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Fig. 2. Top: The map of CMB fluctuations obtained after Wiener filtering of the 
data using parameters estimated with SMICA. Bottom: The map of the second 
component (galactic emission) obtained in the same way. 



